
# coding: utf-8

# # HW1_ZTT
# 

# In[1]:


import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns

get_ipython().run_line_magic('matplotlib', 'inline')

data = pd.read_csv('day.csv')
data.head()


# 这里我把 instant dteday 以及 casual 和 registered 的数据给去掉，因为我觉得特征不相关

# In[2]:


data=data.drop(['instant','dteday','casual','registered'],axis=1)


# In[3]:


data.info()


# In[4]:


data.isnull().sum()


# In[5]:


data.describe()


# In[6]:


# cnt的直方图，感觉看不出什么东西来 - -  这里有个问题怎么判断数据中是否存在负数，图里面出现
#了-2000 是因为有负数的存在吗？ 数据基本和正态分布比较接近
fig = plt.figure()
sns.distplot(data.cnt.values, bins=50, kde=True)
plt.xlabel('Median value of cnt', fontsize=12)
plt.show()


# In[7]:


#cnt的散点图 data.shape 表示行数  这里看着感觉分布比较均匀 需不需要剔除2000以下的值呢
plt.scatter(range(data.shape[0]),data['cnt'].values,color='red')
plt.title('Distribution of cnt')


# 从这个结果上来看，2012年的整体使用量明显要高于2011年

# In[9]:


#季节的直方图,这个好像没什么意义，一年四季一样长啊
sns.countplot(data.season, order=[1,2,3,4]);
plt.xlabel('season');
plt.ylabel('Number of occurrences');


# In[10]:


#月份的直方图，好像也没什么意义
fig = plt.figure()
sns.distplot(data.mnth.values, bins=30, kde=False)
plt.xlabel('proportion of residential land zoned', fontsize=12)
plt.show()


# In[11]:


#天气直方图
sns.countplot(data.weathersit, order=[1,2,3,4]);
plt.xlabel('weathersit');
plt.ylabel('Number of occurrences');


# In[12]:


#气温的直方图 
fig = plt.figure()
sns.distplot(data.temp.values, bins=30, kde=False)
plt.xlabel('proportion of residential land zoned', fontsize=12)
plt.show()


# ## 两两特征之间的相关性

# In[13]:


cols=data.columns

data_corr = data.corr().abs()

data_corr.shape


# In[14]:


#从特征相关性的图来看季节和月份 以及 气温和体感温度 两组的相关性特别高
plt.subplots(figsize=(15,10))
sns.heatmap(data_corr,annot=True)

sns.heatmap(data_corr,mask=data_corr<1,cbar=False)

plt.show()


# In[15]:


#从结果中可以看出 temp 和  atemp 以及 seanson 和mnth 
#的相关性非常高 考虑是否可以剔除atemp 和 season 两组数据
threshold=0.5
corr_list=[]
size=data_corr.shape[0]
for i in range(0,size):
    for j in range(i+1,size):
        if (data_corr.iloc[i,j]>= threshold and data_corr.iloc[i,j]<1) or (data_corr.iloc[i,j] <0 and data_corr.iloc[i,j] <= -threshold):
            corr_list.append([data_corr.iloc[i,j],i,j])
            
s_corr_list = sorted(corr_list,key=lambda x: -abs(x[0]))

for v,i,j in s_corr_list:
    print('%s and %s = %.2f'%(cols[i],cols[j],v))


# In[16]:


data_new=data.drop(['season','atemp'],axis=1)
data_new.head()


# ### 数据准备以及数据分割

# In[31]:


data_train=data_new[data.yr==0]
data_test=data_new[data.yr==1]

y_train=data_train['cnt'].values
y_test=data_test['cnt'].values

X_train=data_train.drop(['cnt','yr'],axis=1)
X_test=data_test.drop(['cnt','yr'],axis=1)

columns = X_test.columns


# In[18]:


X_train.head()


# ### 数据处理

# In[32]:


from sklearn.preprocessing import StandardScaler
#数据标准化

ss_X=StandardScaler()
ss_y = StandardScaler()

X_train =ss_X.fit_transform(X_train)
X_test =ss_X.transform(X_test)

#这里有个问题，对y做标准化之后 那训练好的模型再预测新的一组数据 需要将预测的y也标准化吗？
#还是能再进行反标准化的，这里有点搞不清楚
y_train=ss_y.fit_transform(y_train.reshape(-1,1))
y_test =ss_y.transform(y_test.reshape(-1,1))


# In[20]:


# 数据归一化
# from sklearn import preprocessing  

# X_train = preprocessing.normalize(X_train, norm='l2')

# X_test = preprocessing.normalize(X_test, norm='l2')


# 数据归一化后模型预测的更不好了 所以这里就把它注释掉了

# ### 确定模型类型

# In[33]:


#线性回归

from sklearn.linear_model import LinearRegression

lr = LinearRegression()

lr.fit(X_train,y_train)

y_test_pred_lr = lr.predict(X_test)
y_train_pred_lr = lr.predict(X_train)

#lr.coef_.T 是指的 模型计算后的各个特征权重

fs = pd.DataFrame({'columns':list(columns),'coef':list((lr.coef_.T))})
fs.sort_values(by=['coef'],ascending=False)


# 从特征的权重数值来看，温度 月份 以及风速 气候等影响均较大 其中温度的影响最大

# ### 模型评价

# In[34]:


from sklearn.metrics import r2_score

print('The r2 score of LinearRegression on test is',r2_score(y_test,y_test_pred_lr))

print('The r2 score of LinearRegression on train is',r2_score(y_train,y_train_pred_lr))


# 从R2结果上来看 对train数据的检验还行，但是对12年的结果预测的并不是很好

# In[35]:


f, ax = plt.subplots(figsize=(7, 5)) 
f.tight_layout() 
ax.hist(y_test - y_test_pred_lr,bins=40, label='Residuals Linear', color='b', alpha=.5); 
ax.set_title("Histogram of Residuals") 
ax.legend(loc='best');


# 结果感觉不太好，大致呈正态分布，有点偏右，估计是有噪声值的原因

# In[36]:


plt.figure(figsize=(4, 3))
plt.scatter(y_test, y_test_pred_lr)
plt.plot([-3, 3], [-3, 3], '--k')   
plt.axis('tight')
plt.xlabel('True num')
plt.ylabel('Predicted num')
plt.tight_layout()


# 这个，好像预测的不是很好，初步估计应该是模型复杂程度太高了，需要加一定的惩罚即正则

# In[43]:


#岭回归／L2正则
from sklearn.linear_model import  RidgeCV

alphas = [ 0.01, 0.1, 1, 10,20]

ridge = RidgeCV(alphas=alphas, store_cv_values=True) 

ridge.fit(X_train, y_train)

y_test_pred_ridge = ridge.predict(X_test)
y_train_pred_ridge = ridge.predict(X_train)

print('The r2 score of RidgeCV on test is', r2_score(y_test, y_test_pred_ridge))
print('The r2 score of RidgeCV on train is', r2_score(y_train, y_train_pred_ridge))


# 从结果上来看，岭回归预测的也不是很好，难道是我特征值选的有问题？

# In[44]:


mse_mean = np.mean(ridge.cv_values_, axis = 0)
plt.plot(np.log10(alphas), mse_mean.reshape(len(alphas),1)) 

plt.xlabel('log(alpha)')
plt.ylabel('mse')
plt.show()

print ('alpha is:', ridge.alpha_)


fs = pd.DataFrame({"columns":list(columns), "coef_lr":list((lr.coef_.T)), "coef_ridge":list((ridge.coef_.T))})
fs.sort_values(by=['coef_lr'],ascending=False)


# In[45]:


#岭回归的结果可视化
plt.figure(figsize=(4, 3))
plt.scatter(y_test, y_test_pred_lr)
plt.plot([-3, 3], [-3, 3], '--k')   
plt.axis('tight')
plt.xlabel('True num')
plt.ylabel('Predicted num')
plt.tight_layout()


# In[46]:


# Lasso
from sklearn.linear_model import LassoCV

alphas = [ 0.01, 0.1, 1, 10]

lasso = LassoCV(alphas=alphas)

lasso.fit(X_train, y_train)

y_test_pred_lasso = lasso.predict(X_test)
y_train_pred_lasso = lasso.predict(X_train)

print('The r2 score of LassoCV on test is',r2_score(y_test, y_test_pred_lasso))      
print('The r2 score of LassoCV on train is',r2_score(y_train, y_train_pred_lasso))


# In[47]:


mses = np.mean(lasso.mse_path_, axis = 1)
plt.plot(np.log10(lasso.alphas_), mses) 

plt.xlabel('log(alpha)')
plt.ylabel('mse')
plt.show()    
            
print ('alpha is:', lasso.alpha_)

# 看看各特征的权重系数，系数的绝对值大小可视为该特征的重要性
fs = pd.DataFrame({"columns":list(columns), "coef_lr":list((lr.coef_.T)), "coef_ridge":list((ridge.coef_.T)), "coef_lasso":list((lasso.coef_.T))})
fs.sort_values(by=['coef_lr'],ascending=False)


# 从岭回归和lasso回归的结果上来看，对2012年的预测结果同样不好，结合前面的数据探索，
# 我认为应该是2012年的单车使用人数整体偏高的原因，由此导致用2011的数据训练出的模型
# （无论是线性回归还是岭回归以及lasso回归）来预测2012年的数据，预测值较真实值整体偏低。
